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We report results of a systematic analysis of spatial solitons in the model of ID photonic crystals, 
built as a periodic lattice of waveguiding channels, of width D, separated by empty channels of width 
L — D. The system is characterized by its structural "duty cycle", DC = D/L. In the case of the 
self-defocusing (SDF) intrinsic nonlinearity in the channels, one can predict new effects caused by 
competition between the linear trapping potential and the effective nonlinear repulsive one. Several 
species of solitons are found in the first two finite bandgaps of the SDF model, as well as a family 
of fundamental solitons in the semi-infinite gap of the system with the self-focusing nonlinearity. 
At moderate values of DC (such as 0.50), both fundamental and higher-order solitons populating 
the second bandgap of the SDF model suffer destabilization with the increase of the total power. 
Passing the destabilization point, the solitons assume a flat-top shape, while the shape of unstable 
solitons gets inverted, with local maxima appearing in empty layers. In the model with narrow 
channels (around DC = 0.25), fundamental and higher-order solitons exist only in the first finite 
bandgap, where they are stable, despite the fact that they also feature the inverted shape. 



I. INTRODUCTION AND THE MODEL 



The potential offered by photonic crystals (PhCs) for various applications, as well for the development of funda- 
mental studies is well known [H-[8[. The combination of PhC waveguiding structures with the Kerr or saturable 
nonlinearity may give rise to spatial solitons, that were studied theoretically in various settings [9|-[l2]. Experimen- 
tally, spatial solitons similar to those expected in PhCs have been created in photorefractive media with optically 
induced lattices, using the methods developed in Refs. [l^j-lll]- In addition to spatial solitons, possibilities of the 
creation of temporal solitons were also considered in PhC and related nonlinear media [2^-[23| (the temporal solitons 
in PhCs ma y re semble their counterparts studied theoretically and experimentally in photonic-crystal fibers, see, e.g., 
works. 

The fundamental model of effectively one-dimensional (ID) PhCs is built as a periodic array of material stripes, 
with linear refractive index {nQ)^^^ > 1 and Kerr coefficient (?T-2)inat alternating with empty layers (filled with 
air), that have rig = 1 and n2 = 0. In a more general case, one may consider an "all-solid" model, which assumes 
an alternation of stripes made of different materials. The structure based on the periodic alternation of stripes with 
different local properties is usually called the Kronig-Penney (KP) model [ssj . 

In the case when only no is modulated according to the KP model, while n2 7^ is constant, the model gives rise 
to ID solitons and other nonlinear states that were studied in Refs. [34, _35|. A combination of the linear KP part 
with the uniform cubic-quintic nonlinearity was studied too [3^, js^l ■ A characteristic feature of the latter model is 
bistability of solitons and a plethora of stable multi-peak localized patterns. 

In a realistic PhC model, both the linear and nonlinear local coefficients should be periodically modulated . The 
corresponding nonlinear Schrodinger (NLS) equation for the spatial evolution of the complex amplitude of the elec- 
tromagnetic field, ^'(a;,z), along propagation distance z takes the following normalized form: 

i*. + *.. + I^(2:)(l + a|«'|')«'. (1) 

Here, term "^^^ accounts for the transverse diffraction in the usual paraxial approximation (x is the transverse 
coordinate), and the KP modulation function describes a periodic array of guiding channels, each of width D and 
depth U > 0, which arc separated by buffer (empty) layers of width L ~ D, i.e., D/L = DC may be considered as the 
"duty cycle" of the underlying material pattern: 

'{k"' t <r<^"A;? - ^ 0, ±1, ±2... . (2, 

Coefficient cr = ±1 in Eq. ^ determines the sign of the nonlinearity, a = +1 and —1 corresponding, severally, to 
the self-focusing (SF) and self-defocusing (SDF) material. By means of an obvious rescaling, one can fix L = 27r [and 
also |cr| = 1, leaving the total power, Q = /J*^ |^'(a;)|^ dx, as a control parameter of stationary solutions]. 
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The limit case corresponding to DC ^0, L'^ ^ 00 is the "Dirac comb" , with W{x) replaced by a chain of (5-functions. 
This variant of the PhC model was studied in Refs. [H, Our objective is to consider the regular version of the 
KP system, with an intention to explore effects of the DC parameter on properties of solitons. We report systematic 
results for fixed depth U — 3 (which adequately represents the generic case), and three characteristic values of DC, 
viz., 0.75, 0.50, and 0.25, the second case being the most interestin g on e. The results (the shape of the solitons, 
their stability, etc.) are quite different from those reported in Refs. |l3l . [T3 |. We are chiefly dealing with the SDF 
nonlinearity, a = —1, which is most promising for finding new properties of solitons. Indeed, while the linear part 
of the NLS equation tends to trap them in waveguiding channels, the repulsive nonlinearity produces a competing 
effect, tending to push the field into buffer layers. A noteworthy consequences of the competition is the existence of 
a nontrivial instability border for solitons in the second bandgap, see below. 

In the case of the uniform SDF nonlinearity, which corresponds to Eq. ([1]) with the last term replaced by 
[Vl^(x) — l^'p] the periodic potentia l g ives rise to gap solitons (GSs), which were studied in detail in terms of 
Bose-Einstein condensates (BECs) [slj-jio]. Solitons found in the present work for DC = 0.75 are similar to the ordi- 
nary GSs, while those found for DC = 0.50 and 0.25 are quite different from them. In terms of BEC, it is also possible 
to consider solitons trapped in a purely nonlinear lattice, which corresponds to Eq. |T]) with the last term replaced by 
W{x) l^'l'^ where modulation function W{x) is sinusoidal In the latter context, the competition between 

linear and nonlinear lattices was considered too jH, HJ. However, our results are different from those reported for 
the nonlinear and combined lattice potentials in the BEC models, because in the PhC model the local nonlinearity 
does not change its sign. 

The paper is structured as follows. A brief account of the methods used for the analysis is given in Section II; the 
linear spectrum of the model is also presented in this section. Basic results for the SDF model are reported in Section 
III. These include families of fundamental and multi-peak GSs in the first two finite (alias Bragg-reflection) bandgaps, 
which are always stable in the first bandgap, but may feature an instability threshold in the second. The stability of 
the GSs is identified through the computation of eigenvalues for modes of small perturbations, and is verified by direct 
simulations. Also reported are (weakly unstable) antisymmetric subfundamental solitons (named as per Ref. [45j). 
which are found in the second band gap , and (generally, stable) antisymmetric bound states of fundamental solitons 
(alias "twisted modes" , as per Ref. 22] ) . In Section IV, we present a summary of results for the SF version of the 
model, i.e., the one with cr = +1 in Eq. ([T]). These results are not drastically different from those reported earlier 
in other models combining the periodic potential and attractive Kerr nonlinearity (ssl . [in particular, the stability 
of soliton families found in the semi-infinite (alias total-internal-reflection) gap obeys the known Vakhitov-Kolokolov 
(VK) criterion 47]]. The paper is concluded by Section V. 



II. THE METHODS 



Stationary soliton solutions to Eq. ([T]) are sought for in the ordinary form, '^{x,z) = e''^^$(a:), where A: is a real 
propagation constant, and real function <I>(x) obeys equation 

-k^ + ^o,x + W{x) + a<P^) = 0. (3) 

Equation ([3]) was solved numerically by means of the iterative Newton's method. For finding spatially even solutions, 
the initial guess was taken as $0(2;) = Aosech(aox), with constants ^0 a-nd gq. Odd solutions (the above-mentioned 
subfundamental and twisted modes) were generated by a different initial guess, ^oix) = Aq sin (qox) • sech (aox), with 
another constant qq. Various families of localized solutions arc characterized by the total power as a function of k, 
Q = Q{k). 

The stability of the stationary solutions has been identified through eigenvalues of small-perturbation modes. To 
this end, the perturbed solution was substituted in Eq. ([T]) as {x, z) = e**^^ + u i^) e*'^^ -I- v* {x) e^^ ^] , where 

$0(2;) is a solution of Eq. ([3]), while u{x) and v{x) are complex eigenmodes of the infinitesimal perturbation that 
can grow exponentially at (generally, complex) rate A. The linearization gives rise to the eigenvalue problem in the 
following form: 

-{■&-k)-W{x){l + 2a^l) -aW{x)<^l 

aW ix) $2 _ A:) + {x) (l ^ 2a$2) 

the underlying solution $0 being stable if all the eigenvalues are real. Equations ([4]) have been solved numerically 
with the help of a finite-difference method based on Taylor-series expansions. 

Before proceeding to the presentation of results for soliton families, it is necessary to display the bandgap spectrum 
of the linearized version of Eq. ([31), as the location of solitons should be identified with respect to this spectrum. It has 
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been computed by means of software package SpectrUW (48|. Characteristic examples of the spectrum are displayed 
in Fig. 1 for three different values of the DC. 
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FIG. 1: The bandgap structure, as a function of modulation depth U [see Eq. ([2}], found from the linearization of Eq. ((3)1 for 
(a) D = 3tt/2, (b) D — n, and (c) D = n/2, the respective "duty-cycle" values being DC = 0.75, 0.50, 0.25. Here and in other 
figures that display ranges of k, shaded areas are occupied by Bloch bands. The spectrum was not computed in the black area. 

In the SDF model, soliton families were looked for in the first two finite bandgaps of the spectrum. For very small 
values of the "duty cycle", namely, DC ^ 0.05 and sufficiently large values of U, the results are very similar to those 
reported for the "Dirac-comb" model in Refs. [T3l.[l^. while they are quite different for the above-mentioned values, 
DC ~ 0.75, 0, 50, 0.25 and U ~ 3, for which systematic results are reported below. Naturally, for fixed U and very 
small values of DC, solitons cannot exist with the SDF sign of the nonlinearity. In particular, fixing [/ = 3, we have 
found that the solitons disappear at DC < (DC)^^. « 0.04. On the other hand, in the SF model the solitons do not 
disappear in the same limit, as very narrow solitons with a large total power can be held in a very narrow channel. 

III. THE SELF-DEFOCUSING (SDF) NONLINEARITY 
A. Fundamental and higher-order spatially symmetric gap solitons 

Different families of GSs are identified by the number of peaks in them. In addition to the fundamental (single-peak) 
solitons, we have found families of spatially symmetric (even) localized modes with two, three, and four peaks. Figure 
[2] represents these families by means of respective dependences Q{k), at the three above-mentioned characteristic 
values of DC, in the first two finite bandgaps. In the case of DC = 0.25, the second bandgap is empty, as no soliton 
solutions could be found in it. It is also necessary to mention that, in the latter case, the higher-order solitons are 
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classified as 2-, 3-, and 4-peak ones by analogy with the situation observed for DC = 0.75 and 0.50, while their actual 
shape is different, see Fig. 12] below (in fact, numbers 2, 3 and 4 for these solitons pertain not to the number of peaks 
but rather to the number of channels occupied by each localized mode). 

For DC = 0.75, the GS families are entirely stable, while, in the intermediate case, DC — 0.50, the analysis identifies 
an intrinsic stability border in the middle of the second bandgap. The transition from stable to unstable solutions 
occurs at points where the GSs feature a specific flat-top shape, see below. 




FIG. 2: Total power Q versus propagation constant k for fundamental and higher-order gap-soliton families in the defocusing 
model, is shown at three values of the structural "duty cycle": DC = D/L = 0.75, 0.5, and 0.25, in panels (a), (b), and (c), 
respectively. Families of stable and unstable solutions are shown, severally, by black and red (gray, in the black-and-white 
version) curves. 

For DC = 0.75, typical examples of stable fundamental GSs, together with 2-, 3- and 4-peak solitons, are displayed 
in Fig. [3]and|31 in the first and second bandgaps, respectively. For both DC = 0.75 and DC = 0.50, maxima and 
minima of the local power in stable solitons are always observed, respectively, in waveguiding channels and in empty 
layers between them, while for DC = 0.25 the situation is different, see below. Another relevant conclusion is that 
the A^-peak solitons may be clearly considered as bound states of N fundamental ones, as suggested by the fact that, 
as obvious in Fig. 2, the respective powers are related by Qn ~ -^Qfund- 

In fact, dependences Q{k) and shapes of the solitons in the SDF model with DC — 0.75 are very similar to those 
known for ordinary GSs in the standard model with the KP potential and uniform repulsive nonlinearity, cf. Refs. 
[34, 35, 36]. For DC — 0.5, the fundamental and multi-peak GSs in the first bandgap are close to their counterparts 
displayed in Fig. [3] for DC — 0.75. However, as shown in Fig. O they are essentially different in the second bandgap 
of the model with DC = 0.5 (hence, different from the ordinary GSs too), showing a trend to flattening. 

As said above, a new feature found in the second bandgap of the SDF model with DC = 0.5 is destabilization of 
all branches of the GSs (fundamental and multi-peak ones alike). As seen in Fig. 2(b), this happens at small positive 
values of k. Exactly at the stability border (and close to it, on both sides) all the solitons feature a fiat-top shape, as 
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FIG. 3: Examples of stable fundamental and double-peak solitons (a), and their three- and four-peak counterparts (b), found 
in the first bandgap of the defocusing model with DC = 0.75, for a fixed value of the propagation constant, k — 2. In each panel 
(here and in similar figures below), the right plot shows the spectral plane of stability eigenvalues for the respective soliton, 
A = Ar + iXi (the stability in implied by Ai = 0) . Here and in other figures, the background pattern represents the underlying 
structure of the layered medium. 



shown in Fig. [HI 

Deeper into the instability region, the higher-order GSs develop inverted shapes, with respect to their stable coun- 
terparts, see Fig. [T] on top of the flat-top background, peaks emerge in empty spaces between the waveguiding 
channels. In this case, the fundamental soliton tends to keep a nearly flat-top shape (although its instability is not 
essentially weaker than that of its higher-order counterparts), featuring a shallow depression in the central guiding 
channel. Thus, each localized mode that had N peaks, and — 1 troughs between them, in the stable configuration, 
develops TV — 1 inverted peaks, in locations of the former troughs, after passing the destabilization threshold (which 
corresponds to the flat-top configuration). In fact, this shape transformation suggests that, with the increase of total 
power Q, the effective repulsive nonlinear potential becomes stronger than the trapping linear one, pushing the wave 
into the empty space and destabilizing the trapped state. 
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FIG. 4: The same as in Fig. 4, but in the second bandgap, for k = 1.3. 



When the flat-top fundamental GS becomes unstable, direct simulations demonstrate that growing perturbations 
slowly destroy it, see Fig. 8(a). As for unstable flat-top counterparts of stable multi-peak solitons, the perturbation 
first splits them back into a multi-peak structure, and then gradually destroys them. A typical example of this scenario 
of the instability development is shown, for the GS of the 3-peak type, in Fig. Efb). 

At small values of DC, such as DC — 0.25, GSs of all types are found only in the first bandgap. A peculiarity of 
the soliton families in this case is that they feature peaks in empty layers, being, nevertheless, completely stable, on 
the contrary to their unstable counterparts in the second bandgap found at DC = 0.50, which demonstrate a similar 
structure, as described above. Examples of these stable solitons are displayed in Fig. [9l 



B. Antisymmetric gap solitons 

In addition to the two-peak states displayed above, which may be considered as symmetric bound states of funda- 
mental GSs, a family of their antisymmetric ("twisted") bound states can be readily found too. These bound states 
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FIG. 5: Examples of stable solitons in the second bandgap of the defocusing model with DC — 0.5, for a fixed propagation 
constant, k = 0.4: (a) a fundamental soliton, with Q = 3.20, and a double-peak one, with Q — 6.69; (b) a three- and four-peak 
solitons, with Q — 10.17 and Q = 13.65, respectively. 



turn out to be stable in the first finite bandgap, but unstable in the second, both for DC — 0.75 and 0.50. In the case 
of DC — 0.25, the second bandgap stays empty, as shown above, therefore bound states do not exist in that bandgap 
either, while twisted bound states of fundamental GSs exist and are stable in the first bandgap. These findings are 
presented in Figs. llOl and fTTl An example of the evolution (gradual destruction) of an unstable antisymmetric bound 
state in the second bandgap is shown (for DC = 0.50) in Fig. [T^aV 

It is known that the combination of a periodic potential and spatially uniform SDF nonlinear term (whose power 
may be different from cubic [4^ ) gives rise to an additional family of weakly unstable suhfundamental solitons, in 
the second bandgap only. The solitons were given this name because their total power is smaller than that of the 
fundamental solitons found in the same bandgap |45j]. Each subfundamental soliton is an antisymmetric one, but it 
is not a bound state of fundamental GSs, as it is localized as tightly as fundamental solitons (i.e., it is a twisted 
state squeezed, essentially, into a single guiding channel). Subfundamental solitons can also be found in the second 
bandgap of the present model, for DC = 0.75 and DC = 0.50 (they do not exist at DC = 0.25, as in that case 
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FIG. 6: Weakly unstable nearly flat-top solitons found at A; = in the second bandgap of the defocusing model with CD = 0.50: 
(a) A fundamental soliton, with Q — 5.08, and flat-top counterpart of the double-peak soliton, with Q — 11.33. (b) Flat-top 
counterparts of three- and four-peak solitons, with Q — 17.58 and 23.83, respectively. 



the second bandgap is empty, as shown above). For DC = 0.50, dependence Q(k) for the family of subfundamental 
solitons is shown by the dashed curve in Fig. [TOT a). and an example of the subfundamental soliton is displayed in 
Fig. I12r b). Similar to the ordinary model, the subfundamental solitons in the present model are weakly unstable, 
and spontaneously transform themselves into fundamental GSs belonging to the first bandgap. 



IV. THE SELF-FOCUSING (SF) NONLINEARITY 



In model with SF the sign of the nonlinear term, solitons can be readily found in the semi-infinite gap. As might be 
expected, they are quite similar to solitons previously found in the model with the KP periodic potential and uniform 
nonlinearity In Fig. [131 we display a set of respective curves Q (fc) for different values of DC, and a typical 

example of a stable soliton, for DC = 0.25. The stability of the entire soliton family complies with the prediction of 
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FIG. 7: Unstable solitons with the inverted shape (i.e., with local power maxima in empty layers). All the examples are shown 
for k = —0.22. The top and bottom plots in panel (a) display, respectively, a fundamental soliton, with Q = 9.17 , and the 
former double-peak soliton (which actually features a single-peak structure, as a result of the inversion), with Q = 21.06. The 
top and bottom plots in panel (b) display, respectively, inverted counterparts of former three- and four-peak solitons, with 
powers Q = 30.88 and Q = 42.74. 

the VK criterion [i^l, according to which the positive slope of the curve, dQ/dk > 0, is a necessary stability condition 
(this criterion is irrelevant to SDF models). As concerns the example of the soliton displayed in Fig. [T31 it is worthy 
to note that the soliton remains trapped almost entirely within a single waveguiding channel, even when it is narrow 
(corresponding to DC = 0.25). The same feature has been observed in all other cases. 

In the parameter region investigated in the framework of this work, we have not found multi-peak solitons in the 
semi-infinite gap. Solitons in finite bandgaps of the SF model were not found either (unlike the models combining 
a periodic potential and spatially uniform cubic-quintic nonlinearity, which realizes the competition between SF and 
SDF terms [ai,!!^). 
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FIG. 8; (a) The evolution of an unstable fundamental soliton with k = —0.15 and Q = 17.58 in the self-defocusing model with 
DC = 0.5. (b) The same for an unstable flat-top soliton belonging to the three-peak family, with fc = and Q — 17.58. 

V. CONCLUSION 

The objective of this work was to develop a systematic analysis of spatial solitons in the model of ID photonic 
crystals, built as the periodic alternation of nonlinear waveguiding channels, of width D, and empty layers between 
them, of width L — D. The main parameter of the model is the structural "duty cycle", DC = D/L. We have 
concentrated, chiefly, on the model with the SDF (self-defocusing) sign of the nonlinearity in the channels, as this 
setting makes it possible to study new effects due the competition between the linear trapping potential and its 
effective repulsive nonlinear counterpart. Various types of solitons in this model have been found in the first two finite 
bandgaps (alias Bragg-reflection gaps) of the SDF model. The family of fundamental solitons has been constructed 
too in the semi-infinite gap (alias the one accounted for by the total internal reflection) of the system with the SF 
(sclf-focusing) nonlinearity. In the SDF model with DC = 0.75, and also in the SF one, the solitons are quite similar 
to their counterparts in the respective models with the uniform nonlinearity. On the other hand, the solitons residing 
in the second bandgap of the SDF model with DC = 0.50 strongly differ from the ordinary GSs: both fundamental 
solitons and their spatially symmetric bound states pass a destabilization point, with the increase of the total power. 
At this point, the solitons feature a flat-top shape, while beyond it the shape is inverted, with local maxima emerging 
in empty layers. The shape inversion and destabilization of the solitons are explained by the competition between 
the linear trapping potential and its nonlinear repulsive counterpart. In the system with narrow guiding channels, 
DC — 0.25, GSs are found only in the first finite bandgap, where they are stable (both fundamental ones and their 
bound states), despite the fact that they feature the inverted shape. 

New properties of the fundamental solitons and their bound complexes in the system with moderate or low values of 
DC, 0.50 and 0.25, may find applications to the design of all-optical data-processing schemes based on spatial solitons 
in planar settings. 
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FIG. 9: Examples of stable gap solitons in the first bandgap of the self-defocusing model with DC = 0.25. In panel (a), the 
top and bottom plots represent, respectively, a fundamental soliton (i.e., the one localized in a single channel), with k = —0.1 
and Q = 11.77, and a soliton occupying two channels, with k = —0.05 and Q = 15.26. Panel (b) displays higher-order solitons 
occupying three and four channels, with Q = 23.18 and 31.60, respectively, both pertaining to fc = —0.05. Note that these 
stable solitons feature local maxima of the power in empty channels. 
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